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Abstract 

The construction of stochastic solutions is a powerful method to 
obtain localized solutions in configuration or Fourier space and for 
parallel computation with domain decomposition. Here a stochastic 
solution is obtained for the magnetohydrodynaniics equations. Some 
details are given concerning the numerical implementation of the so- 
lution which is illustrated by an example of generation of long-range 
magnetic fields by a velocity source. 
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1 Introduction 

To solve complex differential problems in large domains, one way to profit 
from parallel computation in multiprocessor machines is to decompose the 
domain into subdomains and assign the problem in each subdomain to a 
different processor. However, one is left with the problem of computation 
of the boundary conditions in the subdomain interfaces. This implies that, 
in addition to the time consumed solving the equation in each subdomain, 
a considerable amount of time will also be consumed in the communication 
between processors. The ideal situation would be to have a method to 
compute local solutions at the interface points without the need for a grid. 
Such a method is implicit in the idea of stochastic solutions. 

In the following the notion of stochastic solution is introduced for linear 
partial differential equations and a method is described that extends this no- 
tion to nonlinear equations. In addition to its use in domain decomposition 
problems, other relevant features of the stochastic solutions are discussed. 
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1.1 The notion of stochastic solution: Linear and nonlinear 
partial differential equations 

Linear elliptic and parabolic equations (both with Cauchy and Dirichlet 
boundary conditions) have a probabilistic interpretation. This is a classical 
result and a standard tool in potential theory. As a simple example consider 
the heat equation 

d 1 5^ 

— ^(t,^) = -^u(t,x) with u{0,x) = f{x) (1) 
The solution may be written either as 

or as 

u{t,x)=E,f{Xt) (3) 

Mx being the expectation value, starting from x, of the Wiener process dXt = 
dWt. 

Whereas Eq.([T]) is a specification of the problem, Eqs.([2]) and ([3]) are 
solutions in the sense that they both provide algorithms for the construction 
of a function satisfying the specification. An important condition for ([2]) and 
([3]) to be considered as solutions is the fact that the algorithm is independent 
of the particular solution, in the first case an integration procedure and in 
the second a solution-independent process. This should be contrasted with 
stochastic processes constructed from a given particular solution, as has 
been done, for example, for the Boltzmann equation. 

Whenever a similar stochastic process algorithm may be associated to 
nonlinear equations, this would provide new exact solutions and new nu- 
merical algorithms. To obtain stochastic solutions for nonlinear equations, 
it is useful to recall that in the linear partial differential equation case the 
stochastic process starts from the point x where the solution is to be com- 
puted and the solution is a functional of the exit values of the process (from 
a space domain D or a space-time domain D x [0,t]). Therefore, it is natu- 
ral to conjecture that for the nonlinear equations the relevant process would 
have a diffusion, propagation or jump component associated to the linear 
part of the equation plus a branching mechanism associated to the nonlin- 
ear part. Then the solution would also be a functional of the exit measures 
generated by the process. 

For the implementation of this conjecture one rewrites the equation as 
an integral one, for which a probabilistic interpretation is given. In the end 
the stochastic solution is the expectation of a functional over a tree-indexed 
measure. The method, which leads to rigorous results, may also be looked 
at, in qualitative terms, as importance sampling evaluation of the Picard 
series. 
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This method was first used in the pioneering paper of McKean [T] for 
the Kolmogorov-Petrovskii-Piskonov (KPP) equation. Later, a similar tech- 
nique was used for the Navier-Stokes [3], the Vlasov-Poisson [1] [5] [6] 
[8], the Euler [7] and a fractional version of the KPP equation [9]. For 
the diffusion equation with (a G [0, 2]) nonlinearities, Dynkin uses a dif- 
ferent method p,0], namely scaling limits leading to superprocesses (for a 
comparison of the McKean- type construction and superprocesses see [llj). 

1.2 Stochastic solutions and numerical algorithms 

Once a stochastic solution is obtained for a partial differential equation, how 
does it stand in comparison with deterministic numerical algorithms? The 
main points to be considered are: 

(a) Stochastic solutions may provide new exact solutions in cases where 
exact solutions were not known before. 

(b) Deterministic algorithms grow exponentially with the dimension d 
of the space, roughly A^'^ (-^ being the linear size of the grid) whereas the 
numerical implementation of a stochastic process only grows with the di- 
mension d. 

(c) Deterministic algorithms aim at obtaining the solution in the whole 
domain. Then, even if an efficient deterministic algorithm exists, the stochas- 
tic algorithm might still be competitive if only localized values of the solution 
are desired. For example by studying only a few high Fourier modes one 
may obtain information on the small scale fluctuations which would require 
a very fine grid in a deterministic algorithm. 

(d) Each sample path of the stochastic process is independent of the 
others. Likewise, paths starting from different points are independent from 
each other. Therefore the stochastic algorithms are a natural choice for 
parallel and distributed computation. 

(e) Stochastic algorithms handle equally well regular and complex bound- 
ary conditions although, of course, the computation of exit times from com- 
plex domains might not be an easy matter. 

(f) Also, as already pointed out, the local nature of the stochastic solu- 
tions make them the most appropriate choice to obtain boundary conditions 
in subdomain interfaces, thus avoiding the time consuming communication 
problem. The remarkable efficiency of this method for domain decomposi- 
tion schemes has been described in |12) |13) |14j . 

2 Stochastic solutions of magnetohydrodynamics 
equations 

Magnetohydrodynamics concerns the dynamics of magnetic fields in electri- 
cally conducting fluids, e.g. plasmas or liquid metals. It is a macroscopic 
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theory which may be considered as an approximation of the Boltzmann's 
equation when the space and time scales are larger than all relevant length 
scales, such as the Debye length or the gyro-radius of the charged particles. 
We consider, in a non-relativistic approximation, the equations of magne- 
tohydrodynamics in 3 dimensions, with non-zero fluid viscosity and electric 
resistivity. The fluid is taken to be incompressible with density p{x, t) = pQ 

constant and uniform. The equations for the velocity V(x,t) of the fluid 

and the magnetic field B{x, t) are: 

= -(V.V)V+ (B-V)B-- VB VP+uV^V+F(x,t) (4) 

ot pofiQ 2/9oAto Po 

_ = -(y.V)B + (B-V)F + — V^S (5) 

ot Po 

V being the kinematic viscosity, r] the resistivity and pQ the vacuum per- 

->■ 

meability. F{x,t) is a forcing term for the fluid velocity. For the Fourier 
transformed quantities, 



v{k,t) = {2tt)-^I'^ j d^xV{x,t)(^^ 
'b{k,t) = {2tt)-^/^ J d^xB{x,t)e''' 
'f{k,t) = {27r)-^/^ J d^xF{x,t)e''' 



use the fact that the divergences of V{x, t), B{x, t), F{x, t) vanish and project 
on the plane orthogonal to the vector k. The projection operator is 

-5> ->■ ^ ^ 

7i"(fc)(0 = ^ - • efc)efc 

with Cjfc = j^. The projection eliminates the gradient terms and, since 
k-v{k,t) = k- b{k, t) = 0, no information is lost on the velocity and magnetic 
fields. Then we obtain 



dv{k,t) 
di 



-uk'^v{k,t) + (27r)-^/2|A;| j d^q \^{q,t) ®k v{k - q,t) 

1 ^ ^ 1 ^ 

b{q,t)®k b{k-q,t)} + (p{k,t) (6) 



PofJ-o 



db{k,t) 



^k^~b{k, t) + (27r)-3/2|A:| / d^q \v{q, t) ®k ^{k - q, t) 



dt Po 

-'b{q,t)^kv{k-q,t)} (7) 
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where (p{k,t) = ttj-^) /(A;, is the Fourier transform of the divergenceless 
part of the forcing and the product ^fc between two vectors ^ , cj is defined 

by 

-S> -5- ~^ -> 

^ w = i{ek ■ O^(fc)^ (8) 

The next step wih be to give a probabihstic interpretation to the mag- 
netohydro dynamics equations by defining a process and an associated func- 
tional that provides the solution. Two constructions will be given. They 
both lead to rigorous results. However, for practical purposes and numerical 
implementation the one that is most convenient will depend on the values 
of the equation physical parameters. 

2.1 Dissipation-controlled stochastic clock 

To guarantee convergence of the functionals associated to the stochastic 

processes it is convenient to rescale the vectors v{k,t) and b{k,t) by a 
/c— dependent function h{k). This rescaling should be familiar from the 
convergence proofs of Picard iteration and will also be called the majorizing 
kernel. Define: 



Xv{k,t) 



v{k,t) 



Xb{k,t) 



b{k,t) 



h{k) ^/pmiHk) 
Then the integral equations equivalent to dMZl) are: 



(9) 



Xv{k,t) 



Xv{k,0)e 



-uk^t 



+ 



+ 00 



ds vk^ e-"'''' I ^'p{k, t-s)+ I d^q 



2 Qv^vvik) Xv{q, t- s)0k Xv{k -q,t- s) + 



2 9v^bb{k)Xb{q, t- s)(S>k Xb{k -q,t- s) 



Hq)h{k - q) 
{h*h){k) 



(10) 



and 

Xb{k,t) 



^ -JLkH , ^ , rj 2 ~^k^s f ,3 h{q)h{k-q) 

Xb{k,0)e ''0 + / ds — k'e / d^q )' \ 

Jo /"o J {h*h){k) 

1 



- gb^vbik) Xv{q, t-s)0k Xb{k - q,t - s)+ 



1 



Ob^bvik) Xb{q, t- s)<^k Xv{k - q,t- s) 



(11) 
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with the functions f?*-!.** and p defined by: 



(k) 

'p{k,t) 



-Qv^bbik) 

-gt^bvik) 

3 ^{k, t) 
vk'^h{k) 



3(27r)-3/2(/i*/i)(fc) 
v\k\h{k) 

2(27r)-3/2^o {h*h){k) 
ri\k\h{k) 



(12) 



The equations (|T0|) . (|TT|) have a probabihstic interpretation as two coupled 
stochastic processes which combine exponential decay and branching. For 



the exponential processes 
abilities up to time t and 



and 



are the survival prob- 



uk' 



ds 



and 



are 



land 



the decay probabilities in the time interval {s,s + ds). On the other 
\^h{q)h{k — q)/{h * h){k) d^q\ is the probability that, given a k mode, one 
obtains a branching to modes q,k — q. The possible branches, three for 
Xv{k,t) and two for Xb{k,t), are selected with equal probabilities (IjIjI) 
and (|, I). The functions g*-^** play the role of coupling constants at the 
branching points and /> is a source term. 

To obtain Xv{k,t) andxb{k,t) the processes are iterated backwards in 
time from time t to time zero. Then, for each realization, starting from 
the values of the initial conditions that are reached at time zero or from 
the source terms, one reconstructs the values at time t following the process 
forward in time and multiplying at each vertex by the appropriate coupling 
constant g^^^^, with the appropriate product (dkik being the Fourier argu- 
ment at that vertex). The solutions Xv{k,t) andxbik,t) are the expectation 
values of this process, obtained by averaging over many realizations. In Fig.l 
we show a typical sample path of the Xbik,t) process. The backwards-in- 
time process, starts from the time t at which the solution is to be computed 
and runs to time zero, except when a source term is sampled, which stops 
that particular branch of the process. 



With the probability structure as defined above, the branching process, 
being identical to a Galton- Watson process, terminates with probability one 
and the number of inputs to the calculation of Xv{k,t) andxbik,t) is finite 
(with probability one). The following bounds are imposed: 

- On the coupling constants: 



g*^**{k) ^ 1 



- On the source term: 



3{2Try^/^{h*h){k) <mm(u,— ) \k\h{k) (13) 



pik,t) 



< 1 



ip{k,t) 



(14) 
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time 



(q,t-s) 



(q-q">s-s") 



(k-q,t-s) 




(q">s-s") 



(k-q-q',t-s-s') 



nr^ — 1 

(q',t-s-s') O 



S{k-q-q",t-s-s"-s"") 



o 6 



Fourier 
mode 



Figure 1: A sample path of the Xb{k,t) process 



On the initial conditions: 



Xv{k,0) 



< 1; 



Xb{k,0) 



< 1 



v{k,0) <h{k); 



b{k,0) 



(15) 

Provided the bounds ()13p - ()15p on the couplings, the forcing and the initial 
conditions are satisfied, the expectation values of Xv and Xb are bounded by 
one in absolute value almost surely. Once a stochastic solution is obtained 
for Xvik,t) and Xb{k,t), one also has, by Q, stochastic solutions for v{k,t) 

and b{k,t). Summarizing: 

Proposition 1 The stochastic processes Xv{k,t) and Xb{k,t), above de- 
scribed, provide stochastic solutions to the magnetohydrodynamics equations, 
existence of such solutions being guaranteed by the bounds dTgjj -dTS )) . 

This being a rigorous result, the mean value over many realizations of 
the process will provide a solution of the magnetohydrodynamics equations. 
However, one notices that if u and — are small, as they indeed are in 
many situations of practical interest, then for short or moderate times most 
process trees have no branching at all. Hence all the contributions related 
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to nonlinear effects are concentrated in just a few exceptional multibranch 
trees. This is a typical large deviation situation, which is even more serious 
than usual because with the < 1 bounds on the couplings we would like to 
be able to neglect the contributions of large multibranch trees. To avoid this 
large deviation problem, when v and ^ are small, the stochastic clock should 
be changed. This we do in the next subsection, where another rigorous 
solution representation is constructed using an externally defined stochastic 
branching clock. 



2.2 Externally defined stochastic clock 

Here we use a time-dependent majorizing kernel 



(16) 



with A a positive real number. The integral equations for these quantities 



are: 



l{k, 0)e-(^+'^'=')* + j^" ds (A + vk") e-(^+'^'=')^ | i^(A;, t - s) + J d^q 
^7^,^^,^,(A;,t - s)^y{q,t - s) (g)k ^{k -q,t- s)+ 



.r^ J i{q)h{k - q) 
{h*h){k) 



■ Iv^bbik, t - s)(ib{q, t- s)i^k %{k - q,t- s) 



(17) 



^ -fb^^bik, t- s) q^{q, t- s)^k 9,{k - q,t - s) + 
^ 7b^bv{k, t - s) %{q, t- s)(g)k ^vik - q,t - s)\ (18) 



with 



3e^^ \k\{h * h)ik) 

"yv^vvykjt) — 'J^^hbykjt) 



76^^,6 (fe,t) = -Jb^bvik^t) 



(27r)3/2 (A + i,k^) h{k) 

2e^^k\{h * h){k) 
(2vr)3/2 (a + ^k^) h{k) 



The stochastic processes associated to (fT7|) and (fTHjl are identical to those 
of (|10p and (jlip . A sufficient condition for convergence of the processes is, 
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as before, assured by keeping the magnitude of all contributions < 1 which, 
together with the fact that the process finishes in finite time with probability 
one, guarantees convergence. This is fulfilled by the following bounds: 
- On the coupling constants: 



3{2TT)-^/'^e^^\k\{h*h){k) < (A + 



mm I jy, -^j k^j h{k) yt 
(20) 



On the source term: 



a{k,t, 



< 1 



ip{k,t) 



- On the initial conditions: 



< 1; 



< 1 



v{k,0) < h{k); 



b{k,0) 



(21) 

< y/Po^o h{k) 
(22) 



In conclusion: 

Proposition 2 The stochastic processes (^y{k,t) and qi)[k,t), above de- 
scribed, provide stochastic solutions to the magnetohydrodynamics equations, 
existence of such solutions being guaranteed by the bounds (2^) -(2 ^) . 

In contrast with the previous solution where the clock is purely dissipation- 
controlled, now the bounds are explicitly time-dependent, meaning that the 
longer t is, the more stringent are the bounds on the majorizing kernel and 
therefore the smaller must be the initial conditions. Therefore this solution 
provides only finite-time solutions and, if longer timer are desired for fixed 
initial conditions, successive finite-time solutions should be patched up. 

This solution provides implementations where the nonlinear effects are 
easier to put in evidence. However, if the dissipation is extremely small the 
solution is not yet fully satisfactory as far as the large deviation problem is 
concerned. This comes about because the most favorable existing kernels 
being those that satisfy either 



or 



{h*h){k) < \k\ h{k) 



{h*h){k) < h{k) 



(23) 



(24) 

inspection of (|2Up implies that it is always min -jj^ that controls the 
magnitude of the kernel. Therefore to handle the extremely small dissipation 
case it would be better to construct a solution that applies also in the non- 
dissipative case {v = r] = Q). Of course, in the non-dissipative case we will 
have the same limitations as in the construction of the finite-time solutions 
of the three-dimensional Euler equation (see, for example Sect. 2. 5 in [15j). 
As before define 



^v{k,t) 



v{k,t) 
e^^h{k) 



%{k,t) 



b{k,t) 



y/PQilQ h{k) 



(25) 
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but we restrict both the initial condition and the kernel to a maximum 
momentum \k\ = Um, that is 



^^ik,0) =<;b{k,0) = , hik) = if\k\>kM (26) 

The integral evolution equations are the same as before ( (jl7p and (jl8|) ). The 
branching probability being 

h{q)h{k-q) 

no new branches are created with \k\ > kM- Therefore this choice of kernel 
effectively projects the equation on the subspace of momentum \k\ < kM- 

Because in the final computation of the functional leading to the solution 

_> 

only the ratios and ^^^^ intervene, consistency is maintained as long 

as the initial conditions satisfy (|26|) . 

In the non-dissipative limit the bound ()20p becomes 



'i{2TT)-'^''^e^\h*h){k)< -r^hik) 

kM 

A general finite-time solution is obtained by a sequence of stochastic pro- 
cesses corresponding to successively larger k'j^ s. 

In the next section we will provide some details on how to implement the 
bounds (fT3|) - p3|) and (p0|) - (p2]) as well as the branching probabilities which 
are then illustrated by an example showing the generation of long-range 
magnetic fields by a velocity source. 

3 Kernel and branching probabilities for the nu- 
merical implementation of the solutions 

3.1 A majorizing kernel and the bounds 

As we have seen before, convergence of the processes requires: 

\k\{h * h){k) < + min {u, jpj h{k) (28) 

(A = for proposition 1). That means 

{h * h){k) < -f' min ( u, — ] \k\h{k) (29) 

for the cases with dissipation and 

\k\{h*h){k) <j"\h{k) (30) 
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in the non-dissipative limit. 

The fohowing majorizing kernel 



satisfies (1291). Indeed, from 



(31) 



3;/ 



1 1 



IT 



one concludes that (1291) holds if 



c < 



7' min ( I/, — 

^ V ' MO 



We now discuss the bounds (fT3|) - p^ for the first solution (proposition 1) 
and (|20| ) - (|22|) for the second solution (proposition 2), using this h (k) kernel. 
For the first case, Eq. (|13p leads to 



(2vr)^ . / rj 
c < ^ — mm z/, — 



(32) 



and from ([12 



\gv^vv{k)\ = \gv^bb{k)\ 
\9b-^vbik)\ = \gb-^bvik)\ 



r,{27r) 



3/2 



(33) 



implying that the coupling constants g are independent of k. Finally, Eqs.([T 
(fT5|) imply 



vik,0) 



- W 



(34) 



b ik,0) 



< 



11 



Choosing the maximum value of c compatible with the bound ()32p 



ISv^** {k) I 
\gb^**{k)\ 

V ik,0) 
t {k,0) 



mm u,— 



< 



< 



< 



„ mm u, — ^ 
3 J- 



K2^)«/^ ■ , 

9^3 mm I u 



' Mo 



ig;^ mm (u, jL 



3;r?F ^^'^I'^'T^ 



(35) 



For the second solution, the bound c' coming from (j20p is 

^3/2 

min ( z/, 



c'<i^min(.,^]e-^* 



(36) 



where t is the maximum time involved in the computation, i.e., the time 
appearing at the l.h.s. of equations (fT7|) . (fTSj) . 

The coupling constants 7 defined in equation p9|) . as well as the source 
term ^, are function of the branching time i — s at which they appear in 
the r.h.s. of equations ([IT]), ([H]). We have then from (fTOD. OD). (|2^ 



ITk— >**(^)^ s)| 
|7f,^**(A;,t - s)| 



(f {k,t — s) 
v{k,0) 



(2^)3/2 A+I^fc^ 



jA(t-s) ^ 



27r'^ 



6 (A;,0) 



< 



< 



< 



(X+uk^) ^X(t-s)j 



(37) 



Choosing as before the the maximum value of c' compatible with the bound 
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([36|1 . we have 



\'yv^**{kit s)| 
\'yb^**{k,t - s)\ 



if {k,t — s) 
v{k,0) 



< 



< 



-Xs 



(38) 



3^ mm ( u. ^ 



Mo 



6 (A;,0) 



^ 3;;3fc2 mm\u,— 



-xt 



In this case the couphngs 7 depend on the momentum of the branching 
particle and on the branching time. Also, one sees that the bound on the ini- 
tial conditions depends on the final time t at which the fields are computed. 
However, when studying the same system for several times, the initial condi- 
tions should be kept fixed. Therefore either one chooses the initial conditions 
to satisfy (|38p for the largest time to be studied or, alternatively, for each 



time a different c' is chosen to satisfy ([37|l . given v {k,0) 
Then, of course, the couplings should be changed accordingly. 

3.2 Branching probabilities 

With the majorizing kernel h (k), the branching probability is 

h{q)h{k-q) \k\ 



and 



b (A;,0) 



p{k,q) 



{h*h){k) 



\1 |2 I |2 3 

\k — q\ \q\ TT-^ 



(39) 



In practice this branching probability is used in the following way: to obtain 
the spherical coordinates of q, namely {\q\,6,ip), pick three independent 

random numbers ri, r2, r3 G [0, 1] and take the direction of k as the reference 
direction. Using the conditional probabilities 



[2 2 

p(0)= p(k,q)\q\ sinOdlql dip = (it — 6) 
J vr^ 

p{q,9) 



p{q\9) 



sm ( 



pie) (7r-e)(fe2 + g2_2|g||A;|cosl 



(40) 
(41) 



and solving 



r2= / p{e')de' 



r3 



P {q,0) dq 
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for 6 and q, one obtains 



2'Kri 

TT (1 



\k\ cos9 + 



r2) 

;| sin tan ((vr 



(42) 



This defines the coordinates relative to A; of a random momentum qq as if k 
were directed along the z— axis. Considering now a matrix that rotates the 
z axis to the direction of fe, e. g. 



D 



one obtains 

For small values of 
mation 



A A 



A 

kx 

yfk'^+k'^ 





kx 

A 

ky 

A 

k. 



(43) 



q = Dqo 



(44) 



(vr — 9) it is convenient to use the following approxi- 



\k\ 



1 - ra 

2]) imply, even for large (vr 



(45) 



This and also the last equation in (j42p imply, even for large (vr — 0), that 
whenever the random number is close to one the momenta of the resulting 
branches are very large. Because large momenta imply short lifetimes of the 
tree branches, one is led to trees with a very large number of branches. 
Computation and generation of such trees is time-consuming. However, 
because the existence bounds p^ - (fT5]) also imply that the contribution 
of each vertex is smaller than one, very large multibranch trees may be 
neglected with a negligible error, unless one is in a large deviation situation, 
as discussed above. If not, an upper bound may be put on the number of 
branches at the stage of tree generation and from the number of neglected 
trees and the upper bound on the number of branches the error is estimated. 



4 Generation of long-range magnetic fields in mag- 
netohydrodynamics 

Here the stochastic solution is illustrated in a situation where one starts 
from a fluid at rest with a very small initial magnetic field and then looks 
for the generation of long-range magnetic fields when the fiuid is driven by 
velocity sources. We take po = iiq = 1, = r] = 0.005. Because of the small 
values of the kinematic viscosity and the resistivity, we will use the solutions 
with externally defined stochastic clock, pT|) and (fT8|) . 
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An important point when studying the generation of long-range mag- 
netic fields by a velocity source is to distinguish the effect of the velocity 
source from the nonlinear transfer of energy between modes. The stochastic 
solution approach is well suited for this study because each evolution tree 
generated by the stochastic algorithm may be computed with sources and 
without sources. Hence, as long as one is looking for the emergence of a 
mode not contained in the initial condition, one may distinguish the effect 
of the source from an eventual nonlinear transfer of energy between modes 
in the absence of the velocity source. This feature would not be so easily 
implemented in other numerical schemes. 

We use the following kernel: 



, _ 2A/27rmin {v, rj/np) _xt 



(A = 0.1) and the initial conditions at time zero 

v{k,0) = 
~b{k,0) = 







if < 0.05 



j:^T^{k)U if A;2>0.05 



with u = y^, "^j' ^ = 0.0001 and 7r(;t) denotes the projection on the 
direction orthogonal to k. 

We force the fluid velocity field at wavenumber a, taking the two follow- 
ing time-independent velocity forcing term, of different helicity: 



F^{r) = 2A (2/3)-^/2 (,os{ay) e~ 



Fy{r) = 2 A (2/3) -3/2 sm{ax) e 4/3 

F^{r) = 2A (2/3)~3/2 [cos(ax) T sin(ay)] e~'^ 
whose Fourier transform is: 



(46) 



fy{k) = -iAe-^^''y+''^' 



/W(fc) = Ae-^'=^{e- 



lb ie 
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The relative helicity Hr of the forcing F is a number between and 1 
which is defined by 



J (frF{r) ■ (V X F){r) 



1/2 

/ d3r|F(r)|2 / d?r\{V x F)(r)|2 



so that Hr{F) = . ^ for the lower sign of the forcing (|46p and Hr{F) = 

for the upper sign. 

We force the fluid at wavenumber a = 5, choosing /3 = 0.185 to have a 
nearly maximum = 0.95 for the lower sign of the forcing. Two different 
source intensities are studied. 

In this setting one then looks for the generation of a magnetic field b{k* , t) 
at k* = (0, 0, 0.1) as the time evolves from t = 0. 

For each generated tree, the forward computations needed to obtain 
the values at time t are performed with the two sources /^^^ and with no 
source. The effects of the velocity source and the nonlinear transfer of energy 
between modes are separated by checking that in some trees the no-source 
result vanishes whereas in others the result is the same with and without 
sources. Of course, there may be cases where one would have both a no- 
source nonlinear effect and a source effect. However we found out that, 
for the range of parameters and times that were used, those situations are 
virtually nonexistent. 

In Figl2] and [3] we show the time evolution of the generated field for 
the two types of sources as well as the contribution of the nonlinear energy 
transfer between modes. Each computed point corresponds to the generation 
of 2 X 10^ trees. The initial conditions are the same for the computations in 
the two figures. The only difference is the intensity of the source. 

One notices that there is a fast growth of the source-generated field for 
small times with a subsequent near saturation for the larger source intensity. 
By contrast the nonlinear energy transfer grows monotonically. One also 
notices that at these parameter values, the helicity of the source does not 
seem to have an observable effect. 

The stochastic process that generates the solution has an essentially mul- 
tiplicative nature. Therefore the reliability of the results should be checked 
by a large deviation analysis rather than by the standard deviation of the 
samples. As shown elsewhere [TB] the large deviation analysis may be done 
by the direct construction of the deviation function from the data. First, 
from the data, one estimates the "free energy" c{t) 

c{t)= lim - log E {exp{tWn)} 
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5 10 15 20 25 



Figure 2: Magnetic field intensity generated by a small intensity velocity 
source with two different helicities (* and A) compared with the nonlinear 
energy transfer (o) 

Wn being the sum of n results. Then one computes the deviation function 

I (x) = sup {tx — c (t)} 
t 

the deviation function giving a logarithmic estimate of the probability P„ of 
a deviation from the sample average 

Pn (dx) X exp {—nl (x)) dx 

For our analysis we have taken 40 samples of 50000 trees each to approximate 
the computation of c{t). The deviation function I (x) is then computed 
numerically. In FigS] and [5] we have plotted the deviation function for a 
sample of 2 x 10^ trees for t = 20 (FigS]) and t = 1 (Figi5]). One notices 
that the results at t = 20 are more reliable that at t = 1, because one sees 
from the values of I{x) that a larger sample would be needed to obtain the 
same kind of relative precision. 

5 Remarks and conclusions 

- Stochastic solutions provide new exact solutions and new numerical algo- 
rithms ^12j [13j [14j . For the particular case of magnetohydrodynamics, the 
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Figure 3: Magnetic field intensity generated by a larger intensity velocity 
source with two different helicities (* and A) compared with the nonlinear 
energy transfer (o) 



local phase-space nature of these solutions may be quite appropriate for the 
studies of plasma turbulence. 

- The convergence bounds for the solutions derived in ()32p to ()38p may 
be quite small. However one should notice that these bounds are too strict 
and obtained in a worst case analysis. In practice, by the very nature of the 
process, the trees finish in finite time (with probability one) and the prob- 
ability of occurrence of many branches in a tree is rather small. Therefore 
much larger values of the parameters may be safely used. 

- In Section 4 we have exhibited the generation of magnetic fields by a 
velocity source and the nature of the process has allowed a clear separation 
of the source effect and the nonlinear transfer of energy between modes. This 
clear separation of effects is a consequence of the nature of the simulation 
method and might be profitably used in other contexts. 

- An important point to take notice of is the need to have an externally 
defined stochastic clock in situations of small dissipation. Otherwise nonlin- 
ear effects are, in practice, very difficult to study in a reliable manner. This 
will also be an important point for the practical application of the stochastic 
solutions for Navier-Stokes [2], [3]. 
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Figure 4: Deviation function for a sample of 2 x 10^ trees at t=20 
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